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Lorentzian signature. Our main result is that the asymptotics for short stochas- 
tic times t encapsulates information also about the equilibrium aspects. The 
moments evaluated with the complex measure and with the real measure de- 
fined by the stochastic Langevin equation have the same t — > asymptotic 
expansion which is shown to be Borel summable. The Borel transform cor- 
rectly reproduces the time dependent moments of the complex measure for all 
t, including their t — > oo equilibrium values. On the other hand the results of a 
direct numerical simulation of the Langevin moments are found to disagree from 
the 'correct' result for t larger than a finite t c . The breakdown time t c increases 
powerlike for decreasing strength of the noise's imaginary part but cannot be 
excluded to be finite for purely real noise. To ascertain the discrepancy we 
also compute the real equilibrium distribution for complex noise explicitly and 
verify that its moments differ from those obtained with the complex measure. 
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1. Introduction 



The complex Langevin method [1, 2] is arguably the best candidate framework to define 
and compute Lorentzian signature functional integrals beyond series expansions. In brief it 
aims at replacing functional averages of some quantity O with the e tS integrand ('complex 
measure') by a limit of averages computed with a real measure on a 'doubled' configuration 
space, see [3] for a review. In the case of a one-component scalar field theory 

where R ti e is a real positive measure for all t and 9 such that R t ((p r ,(pi) — > S(<f> r )S(<f>i) for 
t — > 0. The angle < # < 7r is related to the phase of the Wick rotation and normalized 
such that 9 = and tt correspond to Euclidean and Lorentzian signature, respectively. 
The right hand side can be evaluated numerically for fixed t and 9 by a two component 
version of the usual Langevin method. That is, a pair of stochastic differential equations 
driven by white noise b r , bi is solved for each instance of b r , bi and 0((p b r +i<f>f) is evaluated 
on the solution (j) h r (t) , 4>\{t) at t, after which the ensemble average 0(4> h r yields the 

ratio in (1.1) for fixed t,9. 

The numerical implementation of the Langevin method is seductively simple. The results 
of unguided numerical experiments are however often inconclusive: the large time limit 
of the averages entering the right hand side of (1.1) may fail to exist or may converge to 
the 'wrong' answer, in simple cases where the left hand side of (1.1) can be evaluated by 
other means, see [12] for a recent critical discussion. There is also a considerable body 
of mathematical work on the subject which unfortunately does not seem to cover the 
situations directly relevant to (1.1) and its numerical implementation. What is lacking is 
a theoretical understanding of the domain of validity of the method: if it fails, why docs 
it fail, and conversely in what circumstances can one be assured that the right hand side 
of (1.1) indeed evaluates the left hand side. 

Although the interest in the complex Langevin method comes mostly from field theory 
gradients are not at the core of the issue. In line with earlier investigations [5, 8, 9, 10, 11] 
we will therefore focus on the zero dimensional case and specifically on the paradigmatic 
case of a 4 interaction. To appreciate the origin of the angle 9 and the simplifications 
made compared to a 1+d dimensional 4 theory we quickly run through the main steps of 
the Wick rotation in a lattice formulation. We start from the discretized Minkowski space 
action on a cylinder T x L d , with different lattice spacings ao, a in the temporal and the 
spatial directions. This gives 

Su = yit B A °<^) 2 - l ad ~ 3a ° E( A *^) 2 - t ad ~ la « E *i - r d_la ° £ ^ ' (L2) 

" n n,i n n 
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where n labels the site, T = N a , L = Na, and A o n := n+5 - 0„, A;0 n = - n , 
i = 1,2, Replacing a h-> — ?e~ l6,//2 ao and reverting to lattice units ao = a = 1 the 
partition function becomes 

Z e = / JJ<i0 n exp (So) , 

J n 

So ■= ^/^(A ^^ e «/^[l(AA)^y^ + ^4 (1.3) 

n n 

This multidimensional integral over the real variables (p n is absolutely convergent for all 
— 7r < 9 < 7r. Evaluating the two point function based on (1.3) for A = and 9 = n — e 
with e a positive infinitesimal, one recovers the discretized Feynman propagator. This 
suggests that suitable distributional limits of the Zg based correlators with 9 — > it — e 
would in principle define the Lorentzian signature lattice theory. Setting aside the subtle 
distributional aspects Zg is for 9 ^ also not suited for evaluations based on numerical 
stochastic approaches (Monte Carlo, Langevin, etc) as the exponent is dominated by the 
imaginary parts, inducing an intolerably low signal to noise ratio. The Euclidean signature 
version corresponds to 6 = and circumvents both awkward features, at the expense of a 
more indirect recovery of the Lorentzian signature amplitudes eventually aimed at. The 
gradients in (1.3) are not central to the problem. By discarding them one obtains a zero 
dimensional system whose stochastic quantization resembles 1 + dimensional quantum 
field theory (i.e. quantum mechanics) with the stochastic time providing the added dimen- 
sion. In the following we write S(q) for a (in general complex valued) polynomial action 
for the real variable q which carries a 9 dependence induced by (1.3) interpolating between 
Lorentzian (9 = it) and Euclidean signature (0 = 0). 

Associated with the complex action S(q) is the complex Fokker-Planck equation 

^- Pt ( q ) = P Pt} P = d q (d q + d q S) , (1.4) 

in the stochastic time t with initial condition p t= o(q) = Po(q)- Here P is the transpose of 
the usual Langevin operator L = d 2 —d q Sd q and formally pt(x) = exp(tP)p - One assumes 
that for suitable initial date p t exists and reproduces the complex Boltzmann factor in the 
limit t — > oo, i.e. lim^oo p t = e~ s . For an 'observable' 0(q) the average wrt the complex 
measure p t is defined by 

On the other hand the real measure relevant for the right hand side of (1.1) is defined by 
the following real Fokker-Planck equation 

^R t (x,y)=FR t (x,y), 

P = d X (A R d X - F x ) + dy(A R dy ~Fy), A R ~ Aj = 1 , 

F x = -Re[d x S(x + iy)\ , F y = -lm[d x S(x + iy)\ , (1.6) 



3 



where the forces satisfy d y F x + d x F y = and A R — Aj + 1, Aj > 0, reflects a variant of 
the fluctuation-dissipation theorem. The initial conditions are R (x,y) = po(x)8(y), and 
again one has assume that Rt = exp(£P)i?o is well defined and has a limit for t — > oo. For 
analytic observables (depending on x + iy only) one considers the averages 

/n) JdxdyOjx + iy) R t {x,y) 

{U)Rt fdxdyR t (x,y) ' 1 ' l) 

The Parisi-Klauder conjecture [1, 2] states that under 'suitable subsidiary conditions' both 
averages coincide: 

(0) pt = (0) Rt , for all i>0. (1.8) 

Then lim t ^ 00 ((9)p t = lim t _ > . 00 (C?)jj t should follow and the simulation of the two-component 
Langevin equation can for large t be used to compute the expectation values with the 
complex Boltzmann factor e~ s ^ as in (1.1). 

For definiteness we focus on the case of quartic action S(q) = ae^^q 4 , < a, < 9 < ir, 
where the moments (x p ) pt and {(x + iy) p )n t , with p6N, fully characterize the underlying 
measures. We show in Section 3 that both sets of moments have identical t — > asymptotic 
expansions of the form 

(x P ) Pt ~ E c p , n (-4ae^ 2 )^^~((x + ^n % , c p , n E N . (1.9) 

n>p/2 mod 2 

Moreover the series (1.9) is Borel summable and defines a unique function 

POD 

M p (t) = e-^WU- 1 dse W (-— m )b p {s), < 9 < n , (1.10) 

where 6 p (s) is the Borel sum of (1.9) for a > 0. Third we show 

M p (t) = (x p ) pt , for all t > . (1.11) 

In other words the Borel resummation of the short time asymptotic expansion (1.9) cor- 
rectly captures the dynamics of the complex measure (1.5) aimed at, including its equilib- 
rium aspects. 

Based on (1.9) one might expect that the same holds true for the real measure defined by 
(1.6). However, by direct numerical simulation we find 

M p (t) = ((x + iy) p ) Rt , 0<t<t c {A!), 

M p (t)^((x + iyY) Rt , t>t c (Ai), (1-12) 
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where the 'breakdown' time t c (Aj) depends on the strength A I of the imaginary noise. 
Generically therefore the conjectured equality (1.8) holds for a finite time interval only, 
rendering its use to define the left hand side of (1.1) in terms of the right hand side 
problematic. The breakdown time increases powerlike as A/ ->■ 0. Taking Aj strictly zero 
may lead to an increased sensitivity on initial conditions in the Langevin simulations [6, 9]. 
A refined version of the conjecture (1.8) thus has to read 

(0) pt = lim (0) Rt , < t < lim tciAj) . (1.13) 

For practical purposes (1.13) may suffice provided the temporal variations become small 
before lim j 4 7 _ 5>0 t c (A/) is reached [12]. A theoretical foundation of the method however 
requires a proof that lim^^o t c (Ai) = oo. 

The article is organized as follows. In Section 2 we show that the moments (x p ) pt ad- 
mit a transfer operator representation in terms of a non-selfadjoint propagation kernel 
(e~* H )(g, q') whose properties we examine in detail. In particular (e~* H )(g, q') is shown 
to admit a well-defined spectral representation whose norm-convergence is governed by 
Davies' spectral norms [19, 20, 21]. It follows that the kernel's t — > oo limit is well-defined 
and correctly projects onto its ground state proportional to e -5 ^/ 2 . In Section 3 we derive 
the results (1.9) - (1.11) for the short time asymptotics and its Borel resummation. Section 
4 is devoted to the numerical simulation of the moments leading to (1.12) and the refined 
conjecture (1.13). For Aj > the numerical moments also have finite large t limits which 
however differ from the correct answer lim i _ 5>00 (a; p )p t . To ascertain the disagreement we 
study in Section 5 directly the spectrum and the ground state of the real Fokker-Planck 
operator P. Both are seen to be compatible with the working hypothesis that t i->- e tp 
indeed defines a strongly continuous semigroup with a pointwise non-negative kernel for 
all A j > 0. Its ground state ipo is then used to independently compute the asymptotic 
values lim i _ 5>00 ((a; + iy) p )n t as ((x + iy) p ) vo . Their agreement leaves a non-naive action of 
the semi-group e tpT on holomorphic functions as the likely culprit for the failure of (1.8) 
for Aj > 0. Two appendices contain supplementary material: in Appendix A the noninter- 
acting case is discussed and the complex Langevin method is shown to work perfectly. The 
under-determination of the observable flow mentioned above gives rise to an interesting 
parallelism between the Parisi-Klauder conjecture and quantum mechanical 'supertasks' 
[25] which we describe in Appendix B. 
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2. The transfer operator for the complex sextic oscillator 



The averages based on the complex measure admit a transfer operator realization akin to 
a quantum mechanical system in the stochastic time t: 

/n) _ Jdqdg' Q(q) e^je^, qJe^M) (9 Y] 

{ 1,H jdqdq'e s ^l 2 {e-^)(q,q')e- s ^)^p {q') ' 1 ' ' 

The operators e~ tu , t > 0, generate a semigroup with integral kernel (e~ iH )(g, q') to which 
we will refer to as the 'complex propagation kernel'. Since H is not selfadjoint for 6^0 
none of the usual properties of a transfer operator semigroup can be taken for granted. 
In Section 2.1 we investigate the spectrum of H and in Section 2.2 the spectral norms 
governing its norm-convergence properties for t — > oo. 

We begin by defining H. The complex Fokker-Planck operator P in (1.4) is not symmetric 
even when specialized to a real action. By a similarity transformation it can be mapped 
into a conventional Schrodinger type operator which is selfadjoint for real actions. Starting 
from P in (1.4) we define 

H := - e sM/2 Pe -s( g )/2 ; p e -5( ? ) = ; Ue s( q )/2 = (2 2) 

This gives H = p 2 + Vpp(?) with the "Fokker-Planck" potential 



1\dq) ~2~dq 2 



Vfp = -( — ) (2.3) 



The associated Schrodinger operator factorizes H = (— d q + \d q S){d q + \d q S) and has 
exp{ — \S{q) } as exact ground state with zero energy; a feature related to an underlying 
supersymmetry, see [4, 3] for reviews. 

For selfadjoint Schrodinger operators H = p 2 + V(q), with V(q) a real even polynomial of 
degree 2p, p > 1, the spectrum is known to be positive, purely discrete, and nondegenerate. 
A complete set of real-valued orthonormal eigenfunctions ip n , n > 0, exists and the transfer 
operator (Euclidean signature propagation kernel) has the spectral decomposition 

( e -I H )(g,g') = ^e-^P n (g,g'). (2.4) 

n>0 

The E n are the eigenvalues of H and P n (q, q) = *pn(q)^n(q') is the projector onto the n-th 
eigenspace. In particular the semigroup t >->■ exp{— tH}, t > 0, is strongly continuous and 
the t — > oo limit converges strongly to the projector P onto the ground state. When the 
couplings parameterizing the potential V become complex, the operator H is no longer 
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selfadjoint and none of the above properties can be taken for granted. In fact, in addi- 
tion to the spectrum becoming complex, several other new phenomena occur for complex 
couplings, originally explored by E.B. Davies in the harmonic case [20]. 

In the following we consider the family of hamiltonians 

H = p 2 - e^co 2 q 2 + \e ie q 6 , < cu , < A , < # < vr , (2.5) 
with p = —id/dq. The parameterization is chosen such that the complex quartic action 

S=l>/\e ie ' 2 q A , (2-6) 

gives rise to (2.5) with oo 2 = 3V% and 9 is the Wick rotation angle of the original problem. 
We often keep u 2 as an independent coupling as many structural properties continue to 
hold for generic u 2 . In addition to cu 2 / y/\ = 3 a sequence of other coupling ratios, namely 
uo 2 /\f\ — 4k + 2u + 1, k > 0, v — ±1, gives rise to Schrodinger equations which are 
' quasi- integrable' in the sense that the first k solutions of a given parity v have the form 
ip(q) = Pk-i(q) exp{— \/Ag 4 /4}, with Pk-i a polynomial of degree k—1. Responsible for this 
phenomenon is an underlying dynamical sl2 symmetry, see [16] and the references therein. 
We aim at understanding the structure of the spectrum, of the spectral projections, and of 
the transfer operator associated with the hamiltonians (2.5) as a function of 9. We begin 
with the spectrum. 



2.1 Spectrum 

As mentioned, for 9 = the spectrum of H is positive, purely discrete, and non-degenerate. 
For complex couplings the structure of the discrete spectrum can be understood from a 
scaling argument. Rewriting the eigenvalue equation for 9 = 0, i.e. H]e=o'0 = u > 0, 
in terms of z — exp{— i9/8}q one finds HQ = exp{i9/4}EQ, with Q(z) := ip(exp{i9/8}z). 
Assuming that the angle 9 is restricted such that normalizability is preserved this deter- 
mines the ^-dependence of the eigenvalues. Writing E n (u 2 , A, 0), n > 0, for the eigenvalues 
of H the relation 

E n (tu 2 ,\,9) = e^£„(w 2 ,A,0), to 2 >0, A>0, (2.7) 

links the discrete spectra for real and complex couplings. Normalizability is preserved for 
< 9 < 7r, as anticipated in (2.5): by substitution into the differential equation one sees 
that normalizable wave functions of a \q p potential have a dominant exp{— 2\f\q l+p l 2 / (p+ 
2)} decay. For p = 6 and q = exp{i9/8}z the exponential is damping in z as long as 
j-^g a p OS itive real part. Incidentally (2.7) also provides a means to define the 
discrete spectrum of the 'wrong sign' sextic anharmonic oscillator as \im e ^ n E n (0, A, 9) = 
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^(1 + i)E n (0, A, 0). This leads to a discrete spectrum located on the diagonals of the right 
half plane, in contrast to the E n (0, A, n) = E n (0, —A, 0) which are real but unbounded from 
below. 

To the best of our knowledge the energy levels E n (3V\, A, 0) have not been computed 
before. The factorization H| w= .^ = (— d q + y/\q 3 )(d q + v^Ag 3 ) shows that 

"«M = r<rS^ Al/16exp {-X« 4 }- < 2 ' 8 > 

is the unique normalized ground state with energy Eq = 0. To get the excited state energies 
we use a simple but reliable technique: for a suitable basis on L 2 the matrix elements of 
H define a matrix operator. Diagonalizing truncations of this infinite dimensional matrix 
produces approximate eigenvalues whose accuracy can be tested by probing for truncation 
independence. A natural choice of basis are the Hermite functions. For A > they 
no longer capture the qualitative behavior of the exact eigenfunctions but the resulting 
matrices have a band structure with only a few diagonals populated which is numerically 
advantageous. 

To compute the corresponding matrix elements we express the hamiltonians and the Her- 
mite functions in terms of creation and annihilation operators. With the normalizations 




q = (a* + a)/y/2u, p = i(a* — a)^u/2, u > 0, [a, a*} = 1, the Hermite functions are given 
by \n) = — =a* n |0), a]0) = 0. Converting (2.5) into a normal ordered expression in terms 
of a*, a the matrix elements between the hermite states are readily obtained and read 

(m|H|n> = |(1 - e e / 2 )(2n+l)5 m , n - |(1 + e e / 2 )[y/n(n-l)6 n - 2 ,m + y/m{m-l)8 m . 2 , n \ 

+ (6m-9)5 m _ 4jn + 15(m 2 -m + l)5 m _2,n| + (m -H- n) 

Q (20m 3 + 30m 2 + 40m + 15U m n . (2.9) 

The low lying parts of the spectrum can now be computed by directly diagonalizing the 
hamiltonian matrices truncated to < m, n < N. Stability of the spectrum with increasing 
iV indicates the reliability of the approximative result. 

As a test we recomputed the ground state energy of H with u = 1 for 10~ 3 < A < 10 3 , 
where high accuracy results are available in the literature, see e.g. [15]. For all but very 
large A truncations of iV = 150 are sufficient to obtain the eigenvalues to 6 digits accuracy. 
The results are in perfect agreement with those tabulated in [15]. Note that the widely 
used Hill determinant method occasionally fails for sextic potentials [16]. 

In Table 1 we present results for the low lying eigenvalues in the 'Fokker-Planck' case 
co 2 = 3VX Both the ^-dependence and the A dependence can be extracted analytically, 
see (2.7) and below. The overall structure of the spectrum comes out as 

E n = C n e ie/4 A 1/4 , C n ^C 6 n 3/2 , A>0, 0<#<vr, (2.10) 
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so that it suffices to know the C n . One has C = and for n < 20 a 9 digit accuracy can 
be achieved with truncations iV < 500. We only present the first ten to six digits. The 
limiting behavior as n — > oo follows from a semiclassical analysis, the constant C is known 
analytically, see (2.22). 



n 
C n 


1 

1.935 482 


2 

6.298 496 


3 

11.680 971 


4 

18.042 635 


5 

25.254 605 


n 


6 

33.226111 


7 

41.891010 


8 

51.197908 


9 

61.105 360 


10 

71.579 037 



Table 1: C n for the eigenvalues E n of H with uj 2 = 3\/A- The truncation size is 150 < N < 300, 
the last digits are rounded. 



2.2 The complex propagation kernel 

From the derivation of the phase relation (2.7) one sees that the eigenf unctions fl n of 
(2.5) are related to those, ip n , of the selfadjoint H| e=0 by Q n (q) = t/j n (e ld / 8 q). As a 
consequence the Q n are no longer orthonormal with respect to the L 2 inner product 
(ip,<p) = fdqip(q)*ip(q). Rather the set f2*,f2 n , n > 0, forms a bi-orthogonal basis in 
L 2 [19], 

(n*,a n ) = J d q n n (q)n m ( q ) = s m , n . (2.11) 

The quantities 

(Q n ,Q n ) = Jdqn n (qyn n (q)=:N n (6), (2.12) 
can be interpreted as the norms of projectors P n ip = fi n (fi*,r/>), which satisfy 

PnPrn = S n , m P n , (P^, P m ^) = (if,, Q* n ) (Q n , fi m ) (ST m , ^) , 

(p r ^,p n ^y/ 2 

l<||P„||== S up {M)1/2 =N n . (2.13) 

The inequality follows by specializing to tp = Q n , the Cauchy-Schwarz inequality gives 
||-fn|| < N n , and specialization to ip = fi* enforces equality. The spectral norms (2.12) 
originally introduced by E.B. Davies for the complex harmonic oscillator encode informa- 
tion about the quasi-spectrum and the norm convergence of the heat semigroup generated 
by the non-selfadjoint hamiltonian under consideration. For the complex harmonic oscil- 
lator it was shown in [20] that lim^oo Mn N n = j(9) < oo, with an explicitly known 
constant 7(0). In Appendix A we present a simple generating formula for the N n of the 
complex harmonic oscillator from which the before-mentioned asymptotics can also be un- 
derstood. For anharmonic oscillators it is only known that the N n grow super-polynomially 
[21]. 
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Our goal in the following is to determine the rate of growth of the N n for the Fokker-Planck 
hamiltonian H with oo 2 = 3\/A- In a first step we show that for the class of hamiltonians 
(2.5) the N n only depend on u 2 /V\ and 9, 

2 

N n {u 2 ) \9) = N n {^=,l ) e). (2.14) 

To this end we consider the scaling isometry 

(S x i/,)(q) = A 1 / 1 V(A 1/8 g) , (S x <p,S^) = (</^> , (2.15) 
and note that S\-iqS\ = \~ 1 / 8 q. Hence 

•Sa-iH^a = A 1/4 H , A>0. (2.16) 

/V\, A->1 

Rewriting Hf2„ = E n Q n as S\-iHSx(S\-iQ n ) = E n (S\-iQ n ) and combining (2.16) with 
the fact that the spectrum is nondegenerate gives 

E n (u 2 ,\) = \ 1 ^E n (^=,l), S x -itl n = n n , , n ■ (2.17) 

On the other hand S\ is an isometry, so that 

K,fi„) = l = (^,fi n ) , (2.18) 

which establishes (2.14). 

Next we observe that the pointwise defined spectral sum 

{e-i B )(q,<n:=Y,e-* En Pn(q,<n, (2-19) 

n>0 

with i£ n the eigenvalues in (2.7), is a candidate for the kernel of the transfer operator. 
The iV n 's in principle then are the coefficients in the expansion of the complex partition 
function 

Tr[e _ s H ] 1 ^ e~^ En N n . (2.20) 

n>0 

Since the 'tail' of the partition function corresponding to large quantum numbers n should 
be dominated by semiclassical configurations one expects that the rate of growth of the 
iV n 's can be extracted from a WKB-type evaluation of the partition function. 
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We begin by reconsidering the eigenvalues. For Schrodinger operators p 2 + V(q), with 
V(q) a real even polynomial of degree 2p, p > 1, the spectrum is known to be purely 
discrete and nondegenerate and to scale like E n ~ C n 2p ^ p+1 ^ + ©(n^ -1 )/^" 1 " 1 )), for large 
quantum numbers n, see e.g. [17]. Combined with the scaling law (2.7) one obtains for the 
eigenvalues E n of the hamiltonians (2.5) a scaling behavior 

£ n ~Cn 3/2 +0(n 1/2 ). (2.21) 

In the Fokker-Planck case, u 2 = 3V% it follows from (2.17) that the exact eigenvalues have 
the form anticipated in (2.10), where only the constants C n remain to be determined. One 
expects their large n behavior to be governed by a suitable semi-classical approximation. 
Indeed, application of the SUSY WKB formula [23] gives for the semi-classical eigenvalues 
e n (A) the simple expression 

e n (A) = , C = ( ^ 7 ( / 5 6 / ) 3) ) V4 - 1-14599 . (2.22) 

The A-dependence is evidently of the form mandated by the scaling law (2.17) while the 
n-dependence is in accordance with (2.21). Comparing with (2.10) one sees that SUSY 
WKB yields an expression for the limiting constant lim n _ ) . 00 n~ 3 ^ 2 C n . In fact, the C n 
approach their asymptotic values fairly quickly, rendering (2.22) a good approximation 
to the spectrum. Computing the C ra 's as in Table 1 the ratio C„/(n 3 / 2 C 6 ) comes out as: 
0.854484, 0.999 311, 0.999 993, at n = 1, 10, 100, respectively. The standard WKB formula, 
in contrast, gives the correct n — > oo asymptotics, but a far worse description for small n. 

Motivated by the good description of the spectrum by the SUSY WKB approximation, we 
also consider the associated wave functions. Evaluating the result of [24] in the case at 
hand one finds that the quasi-classical eigenfunctions depend for real A on q only through 
the combination 

z = \ 1/8 n~ 1/4 C~ 1 q , (2.23) 

such that the domain of definition in q corresponds to z G [—1,1] in z. The original norm 
as defined by a g-integral with (n, A)-dependent domain translates into 

2 r 1 27r 1 / 4 n7/6i 3 / 2 

iV " (0) := ^ / 1 < fa!i «M !i »M • T » = rwmNF " lMm ■ (2 - 24) 

The normalizations are such that A^„(0) — > 1 for n — > oo. However the inner products 
with n 7^ m do not approach zero as n + m — > oo with fixed n — m. The quasi-classical 
eigenfunctions for real A therefore do not form an orthonormal set. The explicit expression 
for the n-th quasi-classical eigenfunction comes out as 



n < \ ^ 1/4 A 1/16 Tvr t . . 1 . 

\l n \z) = — , cos — n + nzj(z) + — arcsmz 



(1 - z 6 )V4 



2 J v ' 2 



1 1 7 



f\ z )-=C\F 1 - -, -, -; z b = -0.047460 z b - 0.002 820 z LZ + 0(z ls ) . (2.25) 

ZOO 
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Note that (S\-iQ n )(z) is A-independent, as required by (the quasi-classical counterpart 
of) (2.17). One can now readily restore the ^-dependence and evaluate (2.24) with inte- 
grand |f2 n (,2)| 2 , which defines N n (6) in the quasi-classical approximation. Straightforward 
numerical integration then shows convincingly 

-\nN n (9) ->t(0) < oo. (2.26) 

As an additional test of (2.26) we also investigate the scaling of In N n by direct numerical 
evaluation using an extension of the truncation technique employed for the eigenvalues. 
On account of (2.14) it suffices to evaluate the AT n 's and their asymptotics for one A. Since 
the interplay between the different normalizations is crucial, we spell out the details here: 
inserting a resolution of the identity in terms of real Hermite functions \k), k G No, into 
the eigenvalue equation for H, i.e. 

J2(m\j)Wn) = E n (k\n n ), (2.27) 



one identifies the numerically computed eigenvectors of the truncated hamiltonian matrix 

as 

v^=c n (j\n n ), j = l,...N, neN , c n eC. (2.28) 

Assuming that the exact eigenfunctions Q n form a bi-orthogonal basis normalized according 
to (2.11), (2.12) one expects 

£ka 2 \c n \ 2 N n , 5>W] 2 cl, N^oc. (2.29) 



Numerical diagonalization routines typically produce eigenvectors normalized to have unit 
norm in C N . Based on (2.29) the N n can then be obtained via 



We again consider the Fokker-Planck case oo 2 = 3V\ in detail. From (2.8) one computes 

The N n ,n < no, we compute numerically via (2.30). In a first step we verify the A- 
independence for A^ = 500. Comparing the results for A = 10~ 3 , 10 _1 , 1, 10, 10 3 , one finds 
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Figure 1: lniV n versus linear fit, —1.34 + 0.47n. 



that whenever the N^s are numerically stable they are also A- independent to the same 
accuracy. Specifically, for 9 less than 7r/4 at least N n , n = 0, . . . ,30, are reliable and A- 
independent at the 10~ 5 level, while for larger 9 only N n , n = 0, . . . , 10, or so are reliable 
and A-independent at the same accuracy level. Next we fix A = 1 and evaluate the N n 's 
from N = 1000, 1300, 1500 truncations. The results for 9 = tt/2 are reported in Table 2. 



n 





1 


2 


3 


4 


5 


6 


7 


logiV n 


0.08664 


0.21303 


0.40745 


0.69222 


1.02414 


1.3901 


1.7800 


2.1869 


n 


8 


9 


10 


15 


20 


25 


35 


45 


lniV n 


2.6063 


3.0351 


3.4711 


5.7206 


8.0350 


10.3840 


15.139(1) 


20.105(3) 



Table 2: Logarithm of spectral norms at 6 = tt/2. 



As visible from Fig. 1, a scaling of In N n linear in n is favored already for moderately large 
n, consistent with the previous result (2.26). 

The scaling law (2.26) has implications for the structure of the resolvent and the quasi- 
spectra of H. In our context (2.26) implies that the spectral representation (2.19) is norm 
convergent only if t > t* and divergent otherwise, where using (2.10) 

2 hiiV n . 

** = 777i o max ~rr~ < 00 • ( 2 - 32 ) 

A 1 / 4 COS | n C n 

Importantly the t — > oo limit still projects onto the ground state for all < 9 < n: 

s - lim e -^ {H - Eo) = P , lim e-^-^V = ^oK> V>) > a.e. . (2.33) 

t— too t— >oo 
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This covers the range needed for the Wick rotation and allows one to define the Lorentzian 
signature propagation kernel in terms of 9 — > 7r_. The result is somewhat surprising as 
for 7r/2 < 9 < 7i the potential in (2.5) becomes unbounded from below. In (2.33) the first 
equation implies the second by Cauchy-Schwarz and the strict positivity of |fio(<z)|- In 
extension of the results for the complex harmonic oscillator [20, 21, 22] one expects (2.26) 
also to govern the behavior of (t, 9) >->■ e~* H as a bounded holomorphic semigroup. 

Finally, the analysis leads to an important scaling relation for the averages (2.1). It suffices 
to consider the monomials 0(q) = q p with p even. Writing m p (t) = (q p ) Pt and assuming 
that the initial value distribution po(q) is for complex coupling replaced with po(e ld ^ 8 q) one 
finds 

m p (t) = e-^ 8 m p \ e=0 (e^ 4 t) , (2.34) 
by combining (2.19), (2.16) and (2.7). 



3. Borel resummation of the short time asymptotics 

We now resume the investigation of the putative identity 

(0) Pt = (0) Rt , t>0, (3.1) 

as surveyed in the introduction. The results of Section 2 entail that the left hand side of 
(3.1) is well-defined for all t > and converges to the desired t — > oo average even if the 
underlying quartic action is complex. Our strategy for investigating (3.1) is based on the 
fact that both sides have the same t — > asymptotic expansion, as detailed in Section 3.2. 
For the quartic selfinteraction the coefficients of the expansion turn out to be such that the 
series is amenable to a Borel resummation. The Borel transform defines a unique function 
for all < 9 < tt, which is shown to coincide with the left hand side of (3.1) in Section 3.3. 
Nevertheless this falls short of proving (3.1). In order to disentangle the issues involved 
we recap briefly the heuristic arguments for the validity of the conjecture. 



3.1 Recap of the Parisi-Klauder conjecture 

The rationale for the definition of P in (1.6) can be understood by rewriting it in complex 
coordinates, z = x + iy, z = x — iy. This gives 

P = [d 2 z + d z Sd z + d z (d z S)} + [dl + (d z S)*d, + d- z {d z SY\ + 2(A R + A:)d z d- Z . (3.2) 

We may assume S(z)* = S*(z), where S* has complex conjugated coefficients, e.g. S*(z) = 
^ n>0 a*z n for S(z) = J2 n>0 a n z n ■ One sees the form of P is dictated largely by the 
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requirement that its acts like P on holomorphic functions and maps real functions to real 
functions. The mixing term can formally be removed by taking A R + Aj = 0, i.e. A R = 1/2, 
Aj = —1/2. This however compromises the standard real decomposition of the complex 
Langevin equation and even in the free case leads to equilibrium distributions which are not 
integrable, see Appendix A. In principle the mixed term could be replaced with v(zz)d z d 2 
for any function v of one variable without affecting the reality of the operator or its action 
on holomorphic functions. 

The standard heuristic argument for the validity of (3.1) proceeds by contour deformation: 

p t {x) _ JdyR t (x -iy,y) 



Jdxp t (x) Jdxdy R t (x,y) 

relates both averages provided z R t {z — iy,y) is for fixed j/el analytic with suitable 
fall-off: 



dxdyO(x + iy)R t (x,y) = dy dzO(z)R t (z - iy,y) 

2 JR JR+iy 

= ldxO{x) I dyR t (x-iy,y). (3.4) 
it Jr 

In particular, if R t (x,y) itself has a t — > oo limit ipo(x,y), it should obey 

e~ s{x) Jdy(p (x -iy,y) 



fdxe s ( x ) Jdxdy ip (x,y) 



(3.5) 



For convenience we note the resulting asymptotic form of (3.1) explicitly: there exists an 
integrable nonnegative solution ip Q (x,y) of Fip Q = such that for O of polynomial growth 



(0(x)) e s I (O(z% , (3.6) 



holds, where 



_ JdxO(x)e s{x) _ Jdxdy Q(x + iy) y {x,y) 

Under the assumption that Rt(x,y) has a spectral resolution of the form R t (x,y) = 
J2 n e ~ tEnCni Pn(x,y), F{p n = E n {p n , the large t limit in (1.7) will be dominated by the 
n = term and e~ tE °c will drop out in the ratio (1.7). The mere existence of a limit 
lim t _ > . 00 (C (z))n t = (O(z)) vo therefore only requires ReE n > 0, for all n, not necessarily 
Eq = 0. Likewise the heuristic contour shift (3.4) and its limiting version (3.5) does not 
require E = 0. On general grounds E = must lie in the spectrum of P, see (3.13) below, 
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and E = is also the only ground state energy of P compatible with (3.6). In contrast 
to P in (1.4) where e~ s manifestly is a zero mode, the solutions of Fip = E (p with ip 
integrable and non- negative can in general not be found analytically, for any candidate E , 
so the mere existence of an appropriate tpo in the kernel of P is nontrivial. 

The putative ground state wave function ip has to obey an additional consistency condi- 
tion. Denoting by P T the real adjoint of P one has 

E J dxdy 0(x + iy) ip (x, y) = J dxdy (F t O)(x + iy) Lp (x, y) , 

(¥ T 0)(z) = d 2 z O-d z Sd z O, (3.8) 

using only P<^ = E (f and (3.2) That is: for E ^ the ( ) lf0 averages of a holomorphic 
observable O and its 'dual' Eq 1 F t O have to coincide. For E = the ( ) tpo averages of 
all observables in the image of P T have to vanish; see [12] for an alternative derivation. 

In fact only the Eq = version of (3.8) is compatible with the validity of (3.6). To see 
this, take O(z) = z 2 . Then (3.8) reads (E /2) J dxdy (x+iy) 2 fo(x,y) = J dxdy ip (x,y) — 
f dxdy (zd z S)(z = x+iy) <po(x, y), where by assumption the ( ) (po averages can be replaced 
with ( ) e -s averages. This gives 

^E {x 2 ) e s = 1 - (xd x S) e s . (3.9) 

Upon integrations-by-parts the right hand side vanishes, enforcing E = as the only 
candidate ground state energy for P compatible with (3.6). From a general functional 
analytical principle one can in fact obtain E = irrespective of the validity of (3.6), see 
(3.13). 

Returning to (3.1) and the heuristic argument (3.4) for it, we stress that even in very simple 
interacting theories one has no analytic control over R t (x,y), not even for real arguments 
and the analyticity assumption is little more than a leap of faith. In order to highlight 
the nontrivial nature of the seemingly innocuous steps in (3.4) we spell out some of the 
mathematical underpinnings needed for R t (x,y) to be well-defined. The defining relation 
in (1.6) is a generalized heat equation and the operator P should generate the associated 
semigroup. As such it must satisfy a number of necessary conditions which we first list 
and then comment on: 

(i) F t O(x + iy) = [{d 2 z - d z Sd z )0](x + iy). 

(ii) P generates a strongly continuous semi-group t i->- e' p : L 1 — > L 1 , whose kernel 
(e* p )(x, y\ x', y') is pointwise positive. 

(iii) P has a unique positive ground state ipo G L 1 with zero energy, Fip = 0. 

Condition (i) is necessary for the differentiated version of (3.1) to hold, see Section 3.1. 
It is manifestly satisfied by the proposed real operator P in (1.6), but does not uniquely 
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determine it. Since P is not symmetric the natural functional analytical setting is that of 
a dual pair of Banach spaces, where P acts on one space and P T on its dual. Interpreting 
e tv as an operator on a weighted L l space a sufficient condition for boundedness is 

•ap^/w^-yW^oo, (3,0, 

x',y' J 

for real x,y,x',y' and a suitable weight function w : M, 2 — > K + to be specified later. 
The initial conditions R (x,y) = p (x)5(y) do not lie in L 1 but we assume that R t (x,y) 
is in L 1 and is smooth for all t > 0. The requirement that P is the generator of a 
strongly continuous semigroup poses strong functional analytical constraints, which are 
however indirectly coded in the resolvent and not verifiable by inspection of the differential 
operator. The same holds for the even stronger condition that the evolution kernel is 
pointwise nonnegative. For R t (x,y) both properties in principle follow from the stochastic 
differential equation provided a global solution with the appropriate initial conditions exist. 
A necessary condition for (ii) is that P has a spectrum of the form 

<r(-P) = {E r ± iEi , E r , Ei > 0} , E := inf E r G <r(-P) . (3.11) 

Generally, if P T is the Banach space adjoint of P with respect to a suitable pairing the 
full spectra of P and P T in principle coincide, see e.g. [18], Thm. VI. 7. An eigenvalue E n 
of P however can either be an eigenvalue of P T or lie in P T 's residual spectrum. Here 
we interpret the operators P T and P T as maps from L°° — > L°°, where L°° is dual to a 
weighted L 1 space. For the dual Banach spaces we take 

L 1 ^ 2 , w~ l dxdy) , L°°(M 2 , wdxdy) , w(x, y) = e - 5 ^ 2+y ^ , 5 > , (3.12) 

and similarly for the one variable case. The choice of the weight function ensures that 
all polynomials x p ,y p , (x + iy) p and also e jlX+j ' 2y , e j ^ x+iy ^ for real sources are elements of 
the dual space. The integrability condition on the eigenfunctions of e' p and R t (x,y) is 
stronger than the unweighted one, but should for polynomial actions be easily satisfied. In 
particular the constants are elements of the dual space and are annihilated by P T and P T , 
respectively. Hence G <r(P T ) which explains why zero lies in the spectrum of P. More 
generally one has the following relations among the spectra 

G <r(H) = ct(-P) = a(-P T ) , G a{-F T ) = cr(-P) , (3.13) 

but a{— P T ) C <j(— P T ) does not follow. Note that all spectra refer to different Hilbert or 
Banach spaces. 

Finally, condition (iii) means that the ground state energy E lies in the discrete spectrum 
of the operator. This could be violated in principle but in the present setting of scalar fields 
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with a polynomial interaction it is a feature one would expect to hold. For a symmetric 
semi-group uniqueness of a normalizable ground state follows from strict positivity of the 
kernel in (ii). For a non-symmetric semigroup it hinges on the validity of the t — > oo 
projection property. Since this will be relevant later on let us assume that the kernel in 
(ii) admits a spectral representation of the form 

(e tF )(x,y;x',y') = £ e tE ^ n (x, yty B (s', y') (3.14) 

n>0 

where tp n G L 1 , ip n G L°°, J dxdy <f m (x , y)ip m (x , y) = 5 m>n for m,ti^0, and E = 0. For 
the kernel of e tpT then the same expansion holds with the roles of (f n and ip n interchanged. 
The relevant projection property reads 

Rt(x,y) = J2 etEn ¥n(x,y) dx' Po (x')i> n (x',0) ->• (p (x,y) dx' p (x')i> (x' , 0) . (3.15) 



3.2 Observable flow and its small t expansion 

Under the assumption that Rt{x, y) has sufficient fall-off in x, y one can integrate by parts 
and have the transpose (real adjoint) of P act on the observables. By definition of P its 
transpose acts on holomorphic observables like L z = d\ — d z Sd z . This gives rise to a flow 
equation for the averages (0)n t . Indeed, 

d t JdxdyO{x + iy) R t (x,y) = JdxdyO(x + iy) FRt(x,y) 

= JdxdyF T 0(x + iy)R t (x,y) = J dxdy (L z O)(x + iy) R t (x,y) . (3.16) 

In particular the normalization factor JdxdyRt(x,y) is time independent. The same 
derivation works for averages {0) Pt starting from the complex Fokker-Planck equation 
d t p t = Ppt- Both averages therefore obey itentical systems of 'observable flow equations': 

d t (0} Rt = (d 2 z O} Rt -(d z Sd z O} Rt . 

d t (0} pt = (d 2 x O} Pt - {d x Sd x O) Pt . (3.17) 

The initial data coincide by assumption (0) Ro = (O) po yielding the following simple but 
fruitful Lemma. 

Lemma: The averages (0} Rt and (0} Pt admit asymptotic expansions for t — > which 
coincide 

in\ X^ tn r tn\ r fdx^ n O)(x) Po (x) 
{0) Rt ~^- fn ~{0) pt , c n = jj—^ , (3.18) 



>Rt 

n>0 
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where L = d 2 x — d x Sd x . Generally, to the extent the flow equations (3.17) determine the 
averages they must coincide. 

We add some remarks: (i) The 'observable flow equations' generalize the Schwinger-Dyson 
equation for the partition function in [9, 10]. Only the equilibrium aspects (vanishing time 
derivatives in (3.17)) are utilized in the approach of [9, 10]. (ii) Clearly both averages 
coincide as formal power series in t. In the formal expansion of the R t averages P T acts 
like ~L Z = d z — d z Sd z on O(z) while in the formal expansion of the p t averages P T acts 
like L, both evidently producing the same coefficients. The c n are also the unique solution 
of the recursion relations entailed by inserting a power series ansatz into (3.17). By the 
smoothness property of e ipT and e tpT images for t > the exact averages do admit a series 
expansion in t, and its asymptotic nature readily follows from the uniqueness of the c n . 
(iii) In the noninteracting case detailed in Appendix A the flow equations (3.17) completely 
determine the averages, implying the validity of the Parisi-Klauder conjecture without the 
need to explicitly compute the respective propagation kernels, (iv) In interacting situations 
the flow equation (3.17) will in general fall short of fully determining the flow of the averages 
of generic 0's. Rather upon choosing a 'suitable complete set' of observables the averages 
of a small subset will recursively determine all others. The quartic case detailed later on 
is a good illustration. It does not help that for both averages the time dependence can 
formally be attributed to the observables [12] d t O(t) = (d z — d z Sd z )0{t) differing only in 
the interpretation as functions of z and x; the under-determination persists, (v) A two 
step approach to lift the under-determination is described in [9, 10]. In a first step one 
characterizes solutions of the equilibrium Schwinger-Dyson equations as complexified path 
integrals and then tries to select the one corresponding to the complex Langevin process 
driven by P via minimization of an effective potential. The approach does not discriminate 
between different values of A I > 0. (vi) Our strategy focusses on situations where the 
series (3.21) is Borel summable, i.e. when the c n are sign-alternating and of essentially 
factorial growth. Then both averages can differ only in a fairly prescribed way and the 
Borel transform of the series potentially coincides with the exact result for one or both of 
the averages. 

In the framework of the observable flow the equality (1.8) comes about as follows: 



In contrast to the traditional argument (3.4) the above variant does not require control 
over the the kernel R t (x,y) in the complex plane. The only questionable step now is the 
one marked with ?. 

The transpositions entering the other steps should be unproblematic. For the one vari- 





dxdy [e tpT 0(x + iy)} p (x)5(y) = / dx [e^O^x)} p (x) 





(3.19) 
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able semigroup this is because P is similar to the real symmetric operator H. One has 
_ e s/2 Pe -s/2 = H = H T = -e~ s/2 P T e s/2 , with P T = L x , which implies for the kernels 
(e- m )(x,x') = (e- m )(x',x) and 

(e tF )(x,x') = e s ^ 2 (e- m )(x,x')e~ s ^ 2 = e s ^(e tF )(x', x)e~ s ^ = e tpT (x',x) . (3.20) 

For the two variable case the kernel of the transposed operator is simply (e* pT )(x, x'\ y, y') = 
(e t¥ )(y, y'\ x, x'), and granting the correct domains the second step in (3.21) follows. In 
contrast to P the kernels of e' p and e' p (with the same order of arguments) can however 
not be related by a similarity transformation with a function (multiplication operator). 
Assuming otherwise and inserting an appropriate ansatz one finds Ajd x F y = A^d y F x as 
a necessary condition. Since d x F y = —d y F x follows from the defining relations, this could 
hold only for Ar + Aj = 0. Similarity transformations with other operators are not 
excluded. 

The step ? holds in the sense of an asymptotic series in t on account of the Lemma. For 
the action of the semigroups themselves, with the image supposed to be smooth functions 
in t, we separately highlight the corresponding property: 

(iv) jdx'dy' exp{tF T ){x,0; x' ,y')0{x' + iy') = jdx' exp(tP T ){x,x')0{x'), 
for allt > and all x. 

We regard the invalidity of this extension as the likely culprit for the failure of the method 
whenever it fails. The right hand side is an Aj independent function of x for all t while 
the kernel on the left hand side is Aj dependent. In the t — > oo limit the x dependence 
disappears but a mismatched normalization may have been picked up. The expected 
t — > oo limit of (iv) can be inferred from (3.20) and (3.15) 

0) J dx'dy' ip Q (x\ y'W + iy') = (0} e -s . (3.21) 

Consistency with (3.6) requires that ipo(x,0) = ip > is a constant given by ip Q l = 
J dx'dy' ipo{x', y'). Only if both P and P T have purely discrete spectra does the uniqueness 
of (po imply the uniqueness of ipo and since a constant trivially lies in the kernel of P T 
constancy in y follows as well: ip (x,y) = ip . 

In Section 5 we shall construct ipo explicitly for the quartic selfinteraction and find that 
it depends nontrivially on Aj and so do its averages of holomorphic observables. Our 
failure diagnostics differs from the one in [12], Section 4.2, where the rapid growth of 
/ dx'dy' exp(tF T )(x, y; x', y')0(x'+iy') in y (for the U(l) link model) is argued to invalidate 
the integration by parts in the second equality of (3.19). We essentially define e' p by e' pT 
acting on holomorphic observables so that validity of this step is built in. Table 5 in Section 
5 also provides some direct computational evidence for its legitimacy. Nevertheless (3.21) 
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fails because the ground state ip has the wrong structure. Thus (iv) must fail irrespective 
of any growth property in y or invalid integration by parts. 



3.3 Temporal Borel resummation for quartic actions 

For a polynomial action like S(x) = ax 4 , a = V\e td ^ 2 /2, the powers x p , p = mod 2, are 
a natural complete set of observables. We set 



m p (t) 



real Fokker-Planck evolution , 
complex Fokker-Planck evolution 



(3.22) 



The flow equations (3.17) translate into d t m p = p(p — l)m p _ 2 — 4apm p+2 , p = mod 2. 
Since mo(t) = 1, one sees that all m p (t), p > 4, are determined by m 2 {t) via the recursion 



m r . 



1 

4a 



(p - 3)m p „ 4 



p-2 



d t m p . 



(3.23) 



The first few read 



m±{t) 
m 6 (t) 
m 8 (t) 



8a 



128a 2 
1 



[2 - d t m 2 ] , 



[df m 2 + 96m 2 o;] , 

[dlm 2 + 576ad?m 2 - 960a] 



3072a 3 



(3.24) 



The flow equations (3.17) therefore entail that all moments in (3.22) coincide for all t if the 
second moments coincide, {z 2 )n t = {x 2 ) pt , for all t. This structure gives rise to an amusing 
parallelism to the quantum mechanical supertasks introduced by Norton [25] which we 
outline in Appendix B. 

Either by solving (3.23) or directly from (3.21) one can work out the formal series to any 
desired order. One finds the structure 



m p (t) 



E 



;-4a) 



n-p/2 (2t)' 



n>p/2 mod 2 



c P , n e N . 



(3.25) 



Inserted into (3.23) this gives the recursion relations 

Pip - !; 



c p , p / 2 



-Cp-2,p/2-l 



-ip ~ 3)c p _4,n + — — c p _2,„+i , n>p/2 
p z 



(3.26) 
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In particular c p . p/2 = 2~P/ 2 p\, so that m p (t) = p\t p / 2 /( y p/2)\ + 0(t p / 2+1 ). The coefficients 
relevant for the second and fourth moment are related by c^ n = c 2 , n +i, and the first few 
are listed in Table 3, later on the first nonzero 500 are used. 



C 2 ,l 


1 


C2,3 = C 4 ,2 


6 


C2,5 = Q,4 


216 


C2,7 = C 4 ,6 


22896 


C2,9 — c 4,8 


5360256 


C2,ll = C4,10 


2346299136 



Table 3: Coefficients in the asymptotic series (3.25) for p = 2,4. 

According to the Lemma the series (3.25) should be asymptotic to the exact result for 
t — > 0. By direct summation of partial sums up to 0(t N ) one finds that the result is N 
independent for small < t < U{p), where t*(p) ~ 0.13 for p = 2, 4. For small t the series 
thus defines a function whose values can be compared with the numerical simulations. One 
finds an excellent agreement 

m p (t) | S erics = m p (t) (simulation , < t < t* (p) , P = 2, 4 . (3.27) 

This corroborates the Lemma and more specifically determines an interval < t < (p) in 
which the series (3.25) provides a valid description of the exact functions m p (t) in (3.22). 
Assuming that infpt^Qo) > 0, this also leads to a strategy to prove the Parisi-Klauder 
conjecture for short times based on the series (3.25). 

Next consider the growth rate of the coefficients. With the parameterization 

hg^f = a p (k+ l/2)lnfc-/3 p fc, n = p/2 - 2 + 2k, keN, (3.28) 

71/. 

one finds a p , f3 p close to unity for p = 2, 4 based on the c Pi „ with p/2 < n < 1000. A proof 
for all n can be based on the recursion relations (3.26). Since Inn! ~ (n+ 1/2) Inn — n this 
in combination with the alternating signs indicates that Borel resummation techniques are 
applicable. In choosing the parameterization (3.28) we attributed powerlike terms to a 
redefinition of t via t h-» (Ay/at) 2 . 

With this understanding the conventional Borel sum rather than the Borel-Leroy general- 
ization is applicable. We therefore define 

E n-p/2 (2s) n 
Cp,„(-4a) ^ __ , 

n>p/2 mod 2 

PCX) 

M p (t) := r 1 / dse- s,t b p (s) , (3.29) 
Jo 
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initially for real a > 0. One then finds that the partial sums < n < 1000 produce 
truncation independent results for the Borel sums 62,^4 in the interval < s < 11, as 
shown in the inserts of Figures 2,3. The restricted integration produces a well-defined 
Borel transform in the interval < t < 2, shown in Figures 2,3 (solid lines). The Borel 
transforms extend the low t regime of the original functions in (3.27). and quickly approach 
constant values at around t — 1. Importantly they also agree with the averages defined by 
the complex propagation kernel: 

M p (t) = (x p } pt for all t > . (3.30) 

For a > the equality follows from the uniqueness of the Borel transform subject to (3.28). 
The extension to complex couplings can be done via scaling relations. Generally (3.29) 
implies 

M p (t) = a~ p/4 M p \ a=1 (a 1/2 t), (3.31) 

provided M p (t) is analytic in < argt < tt/4. For the Borel sum analyticity in < arg s < 
7r/4 is manifest and allows one to rewrite (3.31) in a way such that only the Borel sum for 
real arguments enters 

dsexp( y --^ JI )b p (s), (3.32) 

where b p (s) is the Borel sum for a > and M p {t) represents the moments for e ld ^ 2 a. On 
the other hand (3.31) matches precisely the scaling relation (2.34) derived in Section 2 from 
the properties of the complex propagation kernel. This shows (3.30) for all < 9 < it. 

In contrast to the setting in Section 2 the Borel transform (3.32) allows for a direct com- 
putational implementation. The relation (3.32) therefore provides an alternative to the 
complex Langevin method in the case at hand. It directly gives the time dependent mo- 
ments for (2.1) which we will compare in Section 4 with the results obtained from the 
complex Langevin simulations. 

A compelling numerical demonstration of (3.30) for a > is obtained by comparing 
the results of a Id Langevin simulation for (x p ) Pt , a > fixed, with the correspond- 
ing Borel transforms. The results are shown for p = 2,4 in Figures 2,3. One also sees 
that the asymptotic values are approched quickly lim t _ >00 (a; 2 ) pt = {x 2 ) e ~s = 0.47798 and 
lim t ^ 00 (a; 4 )p t = (x 4 ) e -s = 0.5, respectively. 

Finally, the result (3.30) is also consistent with the recursion relation (3.23) remaining 
valid for all t with the m p (t) = (x p ) Pt interpreted as the Borel transform M p (t). Stationary 
distributions have <9fm 2 = 0, k > 1, and (3.23) simplifies to 4am p (oo) = (p — 3)m p _ 4 (oo), 
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Figure 2: Results from Id Langevin simulation at 6 = 0, A = 1 for 7772(f) versus the Borel 
transform M2(t). The insert shows the Borel sum 62 (■?)• 




Figure 3: Results from the Id Langevin simulation at 9 = 0, A = 1 for 777.4(f) versus the Borel 
transform M^(t). The insert shows the Borel sum 64 (s). 

1 

which is readily solved: 

r(2±i) 

m p {oo) = v * y a~ p/4 , p = mod 4, 

Ml) 

m p {oo) = v t| ; m 2 (oo)a (2 ~ p)/4 , p = 2 mod 4 . (3.33) 
Ml) 
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This holds irrespective of the validity of the conjecture (3.1) or its t — > oo limiting form 
(3.6). Assuming that the (0) Pt evolution equation in (3.17) remains valid for all t and 
p the asymptotic values must obey (3.33). The limit 7712(00) remains undetermined in 
agreement with the structure in (3.24). A specific choice for 1712(00) will render (3.33) 
compatible with the asymptotic form of (1.8): the limit \im^ OQ (x p ) pt = (x p ) e ~s is of course 
trivially evaluated and gives 



Numerically, (x 2 ) e -.s = 0.47798 A~ 1/4 , (x 4 ) e -s = 0.5 A~ 1/4 , (x 6 ) e s = 0.7169 A" 3 / 4 . One 
sees that (3.34) matches (3.33) iff m 2 (oo) = oT 1 ^ T(3/4)/r(l/4). Implicitly therefore, the 
Borel resummation described before fixes the parameter 7712(00) undetermined by (3.33) 
to precisely this value. 

In the m p (t) = (z p )n t interpretation of the moments the flow equation for (0)n t in (3.17) 
might likewise be valid for all t but with a m 2 (oo) value different from the one above. If 
so, at least the m p (oo), p = mod 4, and the ratios m p (oo)/m2(oo), p = 2 mod 4, must 
come out as in (3.33). In particular 777,4(00) = 0.3536 \~ x l 2 e~ %e l 2 . 



4. Complex Langevin Simulations: Method and Results 

The complex Langevin process underlying (1.6) is realized as a Wiener evolution in the real 
two dimensional space corresponding to the real and imaginary parts of the integration 
variable. Specifically, for a single integration variable with complexification z = x + iy, 
one generates an ensemble of points in the (x, y) plane by integrating the stochastic pair 
of equations 




a 



p/4 , p = mod 2 , Rea > . 



(3.34) 





(4.1) 



where &r, bj are independent Wiener ("white noise") processes, with 



&fl(*i)&fl(* 2 ) 



2A R 5{t 1 - 1 2 ) , 
2A I 5(t l - 1 2 ) , 



&j(fi)&j(t 2 ) 



(4.2) 
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and A R = Aj + 1 (complex fluctuation-dissipation relation). In practice, we realize this 
process by discretizing the time: 

x(t + At) = x(t)-Re^ dS{X d ^ w) ^At + Ab R , 

y(t + At) = y(t) - Im ^^^ At + Abj, (4.3) 

where Abu (resp. Abj) are Gaussian distributed randoms of variance 2A R At (resp. 2AjAt). 
Two important aspects of our numerical implementation are (see also [13]): 

1. We use continuously (i.e. normally) distributed random variables, rather than dis- 
crete random step variables, to avoid discretization artifacts in the generated en- 
semble. Thus, our ensemble values (starting at the origin x = y = 0, say) are not 
"quantized" but fill out a continuous region in the x-y configuration space. 

2. More importantly, we have found that the Brownian process occasionally wanders 
into regions where the "force" functions F x = Re( dS ^ w ^ ) , F y = Im( 9s ^ ty ^ ) become 
quite large. To avoid losing accuracy in the discretized realization (4.3), we therefore 
readjust the time step At at every update to ensure that the variations Ax, Ay 
remain small. In practice, this is done by choosing 



At = — — (4.4) 

1 + \F X \ + \F y \ 1 ; 

with the nominal time step 5 chosen at some suitably small value (typically, 10~ 5 ). 
If one uses a fixed time step, we have found that it is quite common to obtain 
results which appear to converge, but are simply incorrect as a consequence of rare 
excursions which distort the result due to a loss of accuracy in the discretization of 
the Wiener process. 

We have carried out an exhaustive numerical investigation of the Langevin time dependence 
for both real and complex polynomial actions, specifically for the quadratic ("free") S(z) = 
ue l6 / 2 z 2 , oo > ,0 < 9 < 7r, detailed in Appendix A, and the quartic ("interacting") 
S( z ) = \^\e ie ' 2 z\ A > 0, < 9 < 7r, case. The Langevin time dependence of the second 
and fourth moments O(z) = z 2 ,z 4 was studied by generating a large number (typically 
1-4 million) of independent trajectories in the x-y plane (starting at the origin), using the 
adaptive discretized algorithm described above, and terminating each trajectory when the 
desired Langevin time was reached. This allowed us to study both the short time and long 
time asymptotics of the complex Langevin process. 

For definiteness we present the results for the second moments. The behavior of the fourth 
moments is qualitatively similar but the breakdown times are often signalled by a sudden 
increase in fluctuations size. Figure 4 shows a comparison of the results from the Langevin 
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simulation of the second moments with those from the Borel transform (3.32) for several 
values of Aj and fixed 9 = tt/2. The statistical errors for the Langevin results are smaller 
than the size of the symbols. The putative equilibrium values for large t coincide with 
those obtained by Aarts [14]. 
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Figure 4: Second moments from Langevin simulations (red) for various values of Aj and 9 = tt/2 
compared to the Borel transform (blue). 



One sees that in each case there is good agreement for < t < t c (Aj) but disagreement 
for t > t c (Aj). For the time being we define the 'breakdown time' t c (Aj) informally as the 
smallest time at which both results deviate significantly with respect to the (statistical and 
systematic) errors. This rests on the 'experimental' fact that the transition is relatively 
sharp; a theoretical understanding of this phenomenon is currently lacking. Assuming there 
is a sharp transition t c (Aj) visibly increases with decreasing Aj. For the values considered 
one has t c w 0.16,0.22,0.41,0.67 for Ai = 1,0.5,0.2,0.1, respectively. Tentatively this 
suggests a powerlike scaling 



t c (A/) ~ (Az + a)" 7 , a^0, 7^0.6, (4.5) 
1 
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though other functional forms are compatible with the data as well. The precise form of 
the Ai dependence can of course not be pinned down numerically and neither can be the 
limit \mi Al ^t c {Ai). 

However, the breakdown of the agreement for Aj > shows that at best the refined version 
(1.13) of the conjecture can hold true. A similar conclusion was reached along different 
lines in [12]. In much of the early literature Aj = was used by default. On the other 
hand the heuristic derivations of (1.8) do not discriminate between real and complex noise, 
and neither does the framework of [9, 10]. 

Although the Langevin results for Aj > differ for t > t c (Aj) from the 'desired' results, 
they seem likewise to converge to finite asymptotic values m p (oo)\ Langevin- In line with the 
general setting described in Section 3.1 this should equal the average as computed with a 
positive equilibrium distribution tpo in the kernel of P, see (3.15). In the next section we 
compute f directly from the defining relation P</? = and verify that lim t ^ 00 (z p ) Rt = 
{z p )<p indeed holds for the second moments considered. 



5. Spectrum and ground state of P 

As stressed before the non-naive action of the semi-group generated by the real Fokker- 
Planck operator P on holomorphic functions is the likely culprit for the failure of the 
method. Of course one has very little analytic control over the semigroup in question, so 
one might doubt that P generates a well-defined semigroup at all. In the following we 
present numerical results indicating that at least several necessary conditions for this to 
be the case are satisfied. Along the way we also obtain approximations to the ground state 
wave function of P which allows us to independently compute the relevant moments. 

We continue to treat the anharmonic case as a paradigm, S(z) - 
Fokker-Planck operator P in this case reads 



The real 



P = A R d 2 x + A T d 2 y - F x d x - F y d y - d x F x - d y F y , 



F x = -2y/\ 
F v = -2V\ 



(x 3 — 3xy 2 ) cos ^ + (—3x 2 y + y 3 ) sin ^ 

_ _ _ 

(3x 2 y - y 3 ) cos °- + (x 3 - 3xy 2 ) sin ^ 



(5.1) 



with d x F y + d y F x = 0. The operator has two manifest symmetries 



lrr^A _1 / 8 2:,j/^A _1 / 8 ?/ 



A l/4pj 



A=l 



PI 



x— >— x,y—>— y 



P. 



(5.2) 



which constrain the structure of the spectrum and the eigenfunctions. The first one implies 
that the exact eigenvalues of P scale like A 1 / 4 , just as they do for the complex Fokker-Planck 
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operator P, which in turn follows from (2.17) and the fact that H and P are similar. Since 
P is real for all parameter values the complex conjugate of an eigenfunction will again be 
an eigenfunction with the complex conjugate spectral value. This means the spectrum 
of P lies symmetric to the real axis. It is plausibly purely discrete and the numerical 
results presented below are compatible with the structure (3.11), where \~ l ^E n is A- 
independent and E is non-degenerate. Concerning the eigenf unctions F<p n = E n {p n , the 
scaling relation entails that they can be written as a A-independent function evaluated at 
arguments A 1 / 8 ^, A 1 / 8 ?/, 

<p n (x,y) = <p n \ x=1 (\ l/8 x,\ 1/8 y). (5.3) 

The reflection symmetry ¥\ x ^,^ x ^- y = P entails that the eigenfunctions can be chosen 
to have definite parity, ip n (—x, —y) = 5 n ip n (x,y), 5 n G {±1}- Since we expect the ground 
state to be positive it must be invariant 

<Po(x,y) = <p (-x,-y), (5.4) 

which also ensures that odd observables 0(x + iy) have vanishing expectation value (3.7). 



5.1 Spectrum 

In a first step we now investigate the spectrum of P using a variant of the matrix trunca- 
tion technique from Section 2. Anticipating that zero lies in the spectrum the associated 
eigenvector can be used to approximate ipo and to compute the required averages. To this 
end it is useful to choose real basis functions like the products Hj,{x)Hi{y) of Hermite 
functions in x and y. Reality of the approximate tpo can then be seen directly and only the 
test for positivity requires further analysis. The price to pay is that the matrix elements in 
a real basis are more complicated than in a complex basis adapted to (3.2). We introduce 
a pair of oscillators [a, a*] = 1 = [b, b*] by 

a^ = ^=(a* + a), p x = i^(a* - a) = -id x , 

(5.5) 



y = Tyf^ 6 * + 6 ) ' Py = i \j^^ h * ~ h ^ = ~ id y ' 

with to := v^SA 1 / 4 . Further, we decompose P as follows 



where 



6 6 

-P = A R p 2 x + Aj p 2 y + 2u 2 [{y 2 -x 2 )cos- + 2xy sin -] (5.6) 
to 9 to 9 



X(x, y) = —iAco(x 3 — 3xy 2 )p x , Y(x, y) = —iAco(x 3 — 3xy 2 )p y . (5.7) 
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The symmetry (5.2) now amounts to P^x/^^y/^ = wP| w =i, so that to drops out upon 
insertion of (5.5). One finds in terms of the oscillators 

2 

P= (A i? -2cos-)(2a*a+l) + (4 7 + 2 cos -) (26*6+1) 

uj 2 2 

-(A fi +2cos^)(a* 2 + a 2 ) + (- A 7 +2 cos B -) (6* 2 + 6 2 ) + 4 sin Ua* + a) (6* + 6) 
+1 cos I [X(a, 6) - X(6, a)] + i sin | [F (a, 6) + F(6, a)] , 

where 

X(a, 6) = a* 4 -a 4 +2a* 3 a-2a*a 3 + 6a* 2 + 6a*a+3 - 3(6* 2 +26*6+6 2 + l)(a* 2 -a 2 + l) , 
F(a,6) = (a* 3 + 3a* 2 a+3aV + a 3 + 3a* + 3a)(6*-6) 

-3(a* + a)(6* 3 + 6* 2 6-6*6 2 -6 3 + 36* + 6) . (5.8) 

The basis of Hermite functions in y/oux, y/uuy corresponds to 

\m,n) := —^a* m b* n \0) , a|0)=6|0)=0. (5.9) 
ym!n! 



Finally the matrix elements of (5.8) come out as 

ft. 



2 r 6 

(kl\¥\mn) = (A R - 2 cos -) (2k + 1) + (Aj + 2 cos -)(2Z + 1) 

ou i 2 2 



~{A R + 2 COS -) [y/k(k - 1) 5fc-2,m + - l)5 m _ 2 ,fc] 5«,n 

+ ( - A/ + 2 cos -) [\/Z(Z - 1) <$j_ 2 ,„ + - 1) 5„_ 2 ,z] 5 fc , m 

+4sin-(v / fc<5 fc _i jm + Vm^m-i.t) (y/l$i-i, n + v^n-i,/) 

10 10 

+ g cos - [X M ] + g sin - [7 Himn + Y"« fc , nm ] , (5.10) 
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where 



/ /c! / Tfi\ 

Xkl,mn = $l,n\ f [<5fc,m+4 + 2(k — 2)<5fc i?n+2 ] ~ $l,n\ TT [^m,fc+4 + 2(m — 2)5 mjfc+2 ] 

V '/TT-! V rC! 

+ - l)4,m+2 + (6fc + 3)5fc,m] 

- 3[ V / (/ + 2)(/ + l)^ +2 , n + - l)5i-2,n + (2/ + 1)<J,,„] 

x [A/(m + 2)(m+ l)5 fe ,m+2 - \/m(m - l)5 fc , m _ 2 + 4,m] (5.11a) 
>fci,mn = [VI^-l,n - + lflj+i.n] [a/ (m + 3)(m + 2)(m + l)5 m+3 ,fc 

+ 3v / mTT(m + l)5 m+ i,fc + ?>\fmmb m -x,k + \Jm{m - l)(m - 2)6 m _ 3)jfc ] 

- 3[>/fc(Jfc_i, m + VA; + 14+i,m] [v / (n + 3)(n + 2)(n + l)5,, n+3 



+ Vn + l(n + 3)<J,, n+ i - y/n(n - 2)<^ n _ 1 - y/n(n - l)(n - 2)<J, >n _ 3 ] . (5.11b) 

The matrix elements of the transpose operator P T = A R <9 2 + -Ajc^ + F x d x + i 7 ^, are 
obtained similarily and read 

--(kl\F T \mn) = [A R (2k + 1) + A z (2l + l)]5 km 5 ln 



-A R [^k(k - 1) 5 k -2,m + \/m(m - 1) 5 m - 2 ,k] Si >n 

-A^y/Hl - 1) <J,_ 2) „ + \/n(n - 1) 5 n _ 2i i] <5 fc>m 
16 16 



One can check 



(kl\F T \mn) = (mn\F\kl) , (5.13) 
so for the matrix trunctions the last spectral equality in (3.13) is manifestly satisfied. 

The numerical computation of the spectrum is now straightforward. After conversion of 
the double index k, I — 0, . . . , TV — 1 into a single via % = kN + / + 1, the resulting matrix 
Pij := — (2 / u) (kl\F\mn) , i, j = 1, . . . , N 2 , can be diagonalized numerically. The reliability 
and accuracy of the approximate eigenvalues can be assesed by increasing N. Since some 
of the eigenvalues are complex the ordering is done by modulus, i.e. n > m if \E n \ > \E m \. 
With a C code implementing the Hessenberg transformation running times are 1 — 6 hours 
for N = 80—100 on a conventional laptop with effectively 2G RAM usage. In order to reach 
N = 160 without cumbersome recoding we used a 16G RAM workstation with running 
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times 8 — 140 hours. The approximate spectra stabilize convincingly at the O(10~ 3 ) level, 
from the trend upon varying N between 100 — 160 we guesstimated a systematic error. 
For illustration we report the results for a specific case: Aj — 1, — 0, tt/2 in Table 4. By 
the argument after Eq. (3.12) the lowest eigenvalue must be zero. 



n 


A R = 2,A I = 1,6 = 


^ = 2,^ = 1,0 = ^/2 





0.002(2) 


0.006(6) 


1 


2.0725(2) 


1.921(2) 


2 


5.4920(2) 


5.833(2) 


3 


8.5377(2) 


7.885(2) 


4 


8.314(3) ±i 5.798(3) 


8.590(3) ±i 6.092(3) 


5 


13.962(3) ±i 5.998(3) 


13.125(3) ±i 5.254(3) 


6 


14.310(3) ±i 8.346(3) 


15.074(3) ±i 9.160(3) 



Table 4: Low lying eigenvalues for P with A = 1 based on N = 100 — 160 truncations. 

One sees that the first few eigenvalues are real while higher excited states typically come in 
complex conjugate pairs. Importantly the real parts are nonnegative so that the spectrum 
is compatible with P being the generator of a semigroup with real kernel. The spectrum 
for other parameter values was found to be qualitatively similar and is compatible with 
(3.11). Positivity of the kernel is essential for the interpretation of R t as a probability 
measure but is not manifest from the spectrum. 



5.2 Ground state 

In contrast to the quadratic case the determination of the ground state wave function is 
nontrivial. First note that now exp{— S(x + iy)} is not even an eigenfunction of P. This is 
because although d y acts like id x on holomorphic functions, the replacement is illegitimate 
in the d y F y term. This results in e +5 Pe -5 = -ttmd%S - d y F y = 6y/\e~' w/2 (x - iy) 2 . The 
special real solution (A. 17) generalizes to the interacting case but is again not integrable. 
For generic Ar it is tempting to search for elements in the kernel of the form 

<p (x, y) = exp{-A x* - 3B x 3 y - QC x 2 y 2 - 3D xy 3 - E y 4 } . (5.14) 

By direct computation one sees that no such solution exists for Ar ^ 1/2. Since F(p = 
Eoipo is an elliptic second order differential equation one expects the existence of a solution 
to be determined by general principles. The usual existence theorems however refer to 
bounded domains with Dirichlet boundary conditions. In the situation at hand compact- 
ification of IR 2 to the unit square, say, introduces singularities in the coefficient functions 
of P towards the boundary where the solution is supposed to vanish. In the interior the 
desired solution must be non-negative. Establishing (non-)existence of a solution along 
these lines therefore may be non-trivial. 
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Using the matrix truncation technique the numerical determination of the ground state 
is however straightforward. Once the ground state energy Eq is reliably identified at 
the chosen truncation level N one can compute the associated eigenvector of to high 
accuracy by iterated action of the resolvent. The resulting vector Vi, i = 1,...,N 2 is 
converted into the approximate ground state wave function by the transformation 

N 2 

<fo(x,y) = H k (i)(y/ux)Hui ) (y/uy) , (5.15) 

i=l 

where the cu-dependence is restored and k(i) — (i — 1 — mod(i — 1, N))/N, l{i) = mod(i — 
1, N) converts the indices, with mod(m, n) defined as m modulo n. The required basis of 
Hermite functions is conveniently programmed using the recursion relation. As a test on 
the accuracy of the conversion one can check J dxdy <po(x, y) 2 = u^ 1 ^ vf. 

We computed the eigenvector Vi, i = 1, . . . N 2 , for the same parameters sets as in Table 4 
and converted it into the corresponding ground state wave function <p using (5.15). The 
results are shown in Figure 6. It is numerically nontrivial that they come out strictly 
positive and therefore can serve as a probability measure. 




4 2 -2 -4 4 2 -2 -4 



Figure 5: Ground state wave function for P for truncation N = 140 and parameters Ar = 
2,Aj = 1. Left 9 = 0, right 6 = vr/2. 

For nonzero 9 one sees that the principal axes in Figure 5 is rotated by an angle 0/4. By 
analogy with the non- interacting case described in Appendix B one should not expect, 
however, that ipo differs from (po\e=o only by a rigid rotation. For other values of Aj the 
results for the ground state wave functions are qualitatively similar: a double peak is visible 
whose principal axes is rotated by an angle 0/4. The main effect of decreasing Ai is to 
shrink the extension of the level surfaces in the ^-direction and to raise the summits. This 
is illustated in Figure 6. 
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-3-2-10 1 2 3-3-2-10 1 2 3 

Figure 6: Level surfaces of iV = 140 ground state wave functions for 9 = tt/2. Left Aj = 1, right 
Aj = 0.1. 

Finally we computed the second moments for the same parameter values as in the Langevin 
study. The Langevin moments for large times have previously been obtained by Aarts [14]; 
ours are in full agreement. Here we see that the latter also agree within the errors with 
those obtained from the directly computed putative equilibrium measure <po- 



A 1 


1 


0.5 


0.2 


0.1 


m 2 (cp ) 


0.54-i0.41 


0.48 - i0.30 


0.45 -iO. 23 


0.44 - z0.20 


m 2 (oo) L 


0.56-i0.41 


0.49 - zO.29 


0.45 -iO. 23 


0.44 - z0.20 



Table 5: Asymptotic second moments for = tt/2 as computed from the N = 150 ground state 
(po and from the Langevin dynamics. Both agree within the errors but differ from the 'desired' 
one set by the complex measure: m2(oo) = 0.441596 — iO. 182915. 

Hence 

\im(z p ) Rt = (z p ) V0 , Py? = 0, Vo > , (5.16) 

t— >oo 

holds for the second moments considered. This provides direct computational evidence 
that the interplay between the stochastic differential equations (4.1) and the semigroup 
generated by the real Fokker-Planck operator P is indeed as expected. Further the t — > oo 
limit is compatible with the projection property onto the ground state <po, as stipulated 
in (3.15). Nevertheless (3.6) fails for Aj > which via (3.21) signals that property (iv) in 
Section 3.2 fails. 

1 
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6. Conclusions 



We critically reexamined the complex Langevin method in the paradigmatic case of a 
com pi ex measure. In contrast to earlier studies we focussed on the temporal rather 
than the equilibrium aspects of the stochastic dynamics. Our main result is that the short 
time asymptotics encapsulates information also about the equilibrium aspects. Both sides 
of the conjectured identity (1.8) have identical asymptotic expansions around t — 0. The 
coefficients are such that the series is Borel summable and the Borel transform M p (t) 
of the moments' time series correctly captures the dynamics of the complex measure: 
M p (t) = (x p ) Pt , for all t > 0. The resummation of the t — asymptotic expansion 
therefore provides a practically usable alternative to the complex Langevin method, at 
least for low dimensional systems. 

As to the validity of the conjecture itself, the results of sections 4 and 5 provide a coun- 
terexample for its validity for Aj > 0, in line with [14] and similar results in other sys- 
tems [13]. At best the modified form (1.13) of the conjectured identity can therefore 
hold true. Even for this variant, however, counterexamples can be found. For example, 
for Six) = ^\f\e %e l 2 {x 2 ± oo 2 ) 2 we find for Aj = agreement of the putative equilibrum 
Langevin results with the analytical answer when the potential has a single well (plus sign) 
but not for a double well (minus sign) and 6 ^ 0. See [5] for related numerical results. 
In the oo 2 = case considered here the detailed analysis of the non- self adjoint generators 
H and P (spectrum, spectral expansion, ground state) revealed no pathological features 
of the semigroups e~ m , t > 0, and e tp , t > 0. In particular both semigroups project for 
t — > oo correctly only the respective ground states, e~ s ^ 2 and (po. For e~ <H , t > 0, the 
results of Section 2 come close to a mathematical proof thereof. The deviation from selfad- 
jointness produces nontrivial spectral norms which impact the quasi-spectra and limit the 
extension to a holomorphic semigroup. For t real and large relevant in the present context, 
however, the (e tF )(x, x') convolutions converge to e~ s averages even for couplings where 
the potential is unbounded from below. For the e' p , t > 0, semigroup analytical control 
over the projection property is more difficult, but we regard the computational evidence 
in Section 5 as quite convincing. 

Although both sides of (1.8) have well defined limits, for Aj > the limits are different: 

(x p ) e s = hm (a*> ^ lim (z*>) Rt = (z") Vo . (6.1) 

We isolated as the culprit the limited validity of the identity J dx'dy' exp(tP T )(x, 0; x', y') 
0(x' + iy') = Jdx' exp(tP T )(x, x')0{x') , for all t > and all x. The identity holds 
by construction pointwise in x as an asymptotic series in t and both sides are termwise 
independent of Aj. For t — > oo the identity turns into (3.21) which fails for Aj > 0. It 
does so not because an integrations by part step fails but because the ground state ipo has 
the wrong structure. Since ipo is manifestly Aj dependent the ^/-independence of averages 
of holomorphic observables would be a nontrivial bonus property - which ip simply fails 
to have. 
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Another way to look at the problem is in terms of the under-determination of the observ- 
able flow (3.17). The approach of [9, 10] aims at selecting the solution of the stationary 
equations corresponding to the complex Langevin process driven by P via minimization 
of a suitable effective potential. However the one-parametric family of P's used will in 
general not produce results independent of the parameter. The parameter drops out in the 
short time asymptotics which provides one rationale for the temporal Borel resummation. 
When applicable, the Borel resummation augments the missing piece of information in a 
way compatible with certain analyticity properties in t and the coupling. 

Compared to the field theoretical setting aimed at we ignored mass and kinetic terms. In 
line with earlier investigations [5, 8, 9, 10, 11] we regarded the analysis of a pure potential 
interaction as crucial. The inclusion of subleading mass and kinetic terms should not affect 
the qualitative aspects of the picture obtained. In particular it should be interesting to see 
whether the temporal Borel resummation technique extends to other systems. 
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A: The Langevin method for quadratic complex actions 



In the noninteracting case most of the constituents of the complex Langevin method can 
be computed explicitly, although previously only real noise seems to have been considered 
[7, 6]. Here we allow for Aj > and show that the flow equations (3.17) uniquely de- 
termine the time dependent moments without the need of knowing the real and complex 
propagation kernels. The validity of the Parisi-Klauder conjecture for all Ai > follows. 
The equilibrium distribution of the real measure can be obtained in closed form and the 
numerical evaluation of the equilibrium distribution via the complex Langevin simulation 
shows perfect agreement with the analytical result. For the complex propagation kernel 
we briefly review Davies' construction and present a generating function for the spectral 
norms. 

We begin by describing the relation between the familiar Mehler kernel and the complex 
transfer operator with its nontrivial spectral norms. The hamiltonian is that of the complex 
harmonic oscillator, 

H 9 =p 2 +wVV, u>0, 0<6<n, (A.l) 

where the phase is again normalized so that it corresponds to the Wick rotation angle of 
the original problem. Note that for n/2 < 6 < n the potential is unbounded from below. 
On account of a scaling argument (A.l) has spectrum E n = ue l9 ^ 2 (l + 2n), n G No, with 
eigenfunctions £l n (q) ~ H n {e l6 ^q), where the Hermite functions H n are ortho-normalized 
with respect to the L 2 inner product. The complex eigenfunctions Q n also admit a Fock 
space description in terms of creation and annihilation operators 

Flo 1 6 6 

ae = e t6 V 4 . / q _[_ i e ~ %e ^ p = cos -a + i sin -a* , 

V2 y v 7 ^ 2 2 

16 6 

— ie~ %e ^ . — p = cos -a* + i sin -a , (A.2) 
V2^ 2 2 K ' 

where [ae, ae] = [a, a*] = 1. Since a* e = a-e the L 2 adjoint does not preserve the commuta- 
tion relations. Nevertheless He = ue ld (2aeae + 1) holds and \n) := n\~ 1 l 2 a e x |0), ao\0) = 0, 
are the realization of the eigenfunctions fl n . Note that the (^-dependent) |0) can be viewed 
as a coherent state over the 6 = Fock vacuum |0). To define the relevant inner product 
we introduce a linear anti-involution i on the algebra A generated by the ag, de by 

t(ae) := de , i{zx) = zi(x) , i(xy) = i{y)i{x) , z e C, x,y G A . (A. 3) 

Then {x\0),y\0)) 6 := (\0),i(x)y\0)) e with (|0), |0)) e = 1 defines a R-bilinear positive 
quadratic form over A. In particular the Fock space realization of the Q n is an orthonormal 
basis with respect to ( , )e- Returning to position space and the usual L 2 inner product 
{ip, tp) — J dqip(q)*(p(q), this means the fl n form a bi-orthogonal basis in L 2 while the 
quantities {Q n ,Q n ) = N n (9) are the spectral norms, see Eqs. (2.11) - (2.13). On account 
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of (A. 2) they can in principle be computed algebraically. More conveniently a generating 
functional for the N n can be obtained from the Mehler kernel. In the present conventions 
it reads 

V ' V27rsinht I 2tanht sinhti 

= J2^ t{1/2+n) H n (q)H n (q'), (A.4) 

n>0 

where we take u — 1 from now on. Pointwise the substitution q i— >■ e ld ^q, q' \- >■ e l6 ^q' is 
legitimate producing the transfer operator of the complex harmonic oscillator 



(e-§ H °)(e^,e^y) = £ e -^ 2 +")p n (g, g') = exp { -V^ 2 H,}(g, g') , (A.5) 

n>0 



with P„(g, g') = f2 n (g)f2 n (g'). Assuming that for q = q' integration and summation can be 
exchanged one obtains 

/ \ ^ 2 

V(r + r 1 )cos|-2; ^ 



In particular 



1 1 5 — cos# . , . 

cos ' ^ cos^ £ 4cos t,// ^ 



2 2 2 



from which one also infers a powerlike divergence as 6 — > 7r_. In addition one has the 
asymptotics 

lim -lnJV n = 7(0) < oo , (A.8) 

rt— >oo 

where an explicit expression for 7(0) is known [20]. Norm convergence of the sum in (A.5) 
requires that |ReP n — In N n is negative for large n, leading to the conclusion that 

t > U = , (A.9) 

COS 2 

suffices for norm convergence. In particular the t — > 00 limit still projects onto the ground 
state, see (2.33). The semigroup t — > e m in fact defines a bounded holomorphic semigroup 
for certain ^-dependent sectors of the complex t plane [21, 22]. In addition (A. 3) affects 
resolvent estimates and the quasi-spectra [20, 22]. 
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The hamiltonian (A.l) is up to an additive constant the Fokker-Planck hamiltonian for 
the complex action S(q) = Lue t9 ^ 2 q 2 . By construction it is up to a sign isospectral to the 
complex Langevin operator P, viz 

P = d 2 x + 2ue w ' 2 xd x + 2coe w ^ 2 , 

H = -e s ^/ 2 Pe- s ^ 2 = -d 2 x + u 2 e* d x 2 -u, (A.10) 
where we omit the subscript 9 from now on. The real Langevin operator reads 

9 . 9\„ _ / 9 . 0\„ 



P 



A R d 2 + Aid 2 + 2u (x cos - — y sin - ]d x + 2o>( y cos - + a; sin - j d y + 4a; cos - . (A. 11) 

V jLj Z / V Z Z / _ 



The propagation kernel t h- >■ e* p has been studied before in terms of H. The propagation 
kernel t h- >■ e tp is nontrivial even in the free case. Remarkably neither of these objects is 
needed to prove the Parisi-Klauder conjecture for the quadratic action. 

The flow equations (3.17) for the powers (3.19) now assume the simple form d t m p = 
p(p — l)m p _ 2 — 2upm p , for 9 = 0, and the dependence on the phase can trivially be 
restored by rescaling u h-> e %e l 2 uj. These are decoupled differential equations with solution 

m p (t) = e- 2p ^m p {0) + p(p -1)1 dt ' m v _ 2 (i!)e 2 ^ ■ "*> . (A.12) 

Jo 

Since m (t) = 1, all moments are determined recursively up to their values at t — 0. 
Since by assumption the initial values in both interpretations (3.19) coincide, all moments 
coincide for all t, thereby verifying the Parisi-Klauder conjecture in this case. Alternatively, 
one can use the evolution d t O = (d 2 — d z Od z )0 to define 'eigen-observables' and reach the 
same conclusion. The eigenvalue equation translates into 

(d 2 x - 2uxd x )O n = -2nuO n , (A.13) 

which is the defining relation for the n-th Hermite polynomial related by O n = VL n e~ s l 2 to 
the Hermite functions Q n . The flow equations (3.17) are trivially solved 

(O n ) Rt = e- 2 ™\O n ) Ro = e- 2 ^(O n ) P0 = (O n ) pt , (A.14) 

and are equivalent to (A.12). 

Since (A.12) verifies the Parisi-Klauder conjecture one expects that its t — > oo limit verifies 
the asymptotic form of the conjecture. Indeed, for t — > oo the dependence on the initial 
values drops out in (A.12) and the limits are directly related by the recursion 2o;m p (oo) = 
(j) — l)m p _ 2 (oo). Hence 

p\ / 1 \p/2 
m p (oo ' ' 



(P/2) 



^y(^) P =(x p } e -s = (z% , (A.15) 
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where the last equality follows from the uniqueness of the solution (A. 12), assuming only 
that exists. 

It is instructive to see that the equilibrium distribution ip can in the harmonic case be 
found explicitly. We begin by searching for solutions of Py? = E (p in the form 

<p (x,y) = exp{-Ax 2 - 2Bxy - Cy 2 } . (A. 16) 

For generic Ar = Aj + 1 one finds three solutions. In addition to the expected 

<p (x,y) = exp-S(x + iy), E = 2ue- ld/2 , (A.17) 

and its complex conjugate, there is a real solution with E = 0: 

An 



2u cos | (2 cosh § - 1 - cos 9) 
cosh a — cos 9 

w(sin| + sinf ) 
cosh a — cos 9 
2u cos |(2 cosh f + 1 + cos 9) 



Co = iV r 2 ' (A.18) 

cosh a — cos 9 

where we set A R = cosh 2 a/4, A t = sinh 2 a/4. For a = this reduces to the result in [7, 6]. 
The solution becomes non-normalizable for 9 — > 7r_ as <p (x,y) — > 1. After computing the 
averages, however, the limit (9 — > 7r_ is well-defined, see (A. 22) below. The coefficients in 
(A.18) are such that 

Ao Bq \ ( cos - sin - \ ( A + \ / cos j — sin : 



Go r'„ / V -sinf cos f J V A_ J V sinf cos f 



2wcos| 111 . . 

A ± = 2—^, — = -. (A.19) 

cosh | ± cos | A + A_ u; 

Note that the eigenvalues X± are ^-dependent, so the < 9 < it distribution differs from 
the one at 9 = not just by a rigid rotation. 

Specifically for A R = 1/2, i.e. a = in mod i2ir, a second real solution with E = exists: 

<p (x,y) = exp{-S(x + iy)-S*(x-iy)}, Ar — 1/2 . (A.20) 

The solutions (A.17), (A.20) are not in L X (1R 2 ) while (A.18) is integrable for all < 9 < it, 
provided 

9^aeR, or Ann - 9 < ia < Ann + 9 , n E N . (A.21) 



40 



Likewise the constraint (3.5) is satisfied for (A. 17) only with a divergent constant of pro- 
portionality, while for (A. 18) it is satisfied with a finite one. The special solution (A. 20) 
satisfies (3.5) with a finite constant of proportionality only if formally n < 9 < 2ir is 
assumed. This leaves only (A. 18) as an acceptable L 1 (M 2 ) solution. The average {O) vo 
computed from (A. 16), (A. 18) with O(z) = expjz correctly evaluates to 

id 12 

(e^ = exp{^-j 2 } = (e^ , (A.22) 

where Qq{x) = e~ s ( x > is the ground state of P. As a final test we also computed the equi- 
librium distribution of the Langevin dynamics numerically and found excellent agreement 
with the analytical formula for (p Q , see Figure 7. 




Figure 7: Ground state wave function for P with parameters Aft = 2, Aj = 1, oj = 1, 9 = ir/2. 
Left: contour plot of analytical ipo from (A. 16), (A. 18). Right: results from complex Langevin 
simulations. 
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B: Resummation as a quantum mechanical 'supertask' 



The recursion relation (3.23) for the moments is virtually identical to a system of differential 
equations devised by J. Norton in the context of quantum mechanical 'supertasks' [25]. 
This gives rise to an amusing parallelism which we cannot resist mentioning. In brief, 
Norton considers a Schrodinger equation id t ip = Hip, where ip(t) is expanded with respect 
to some L 2 orthonormal basis 

^W = ^/n(t)e i Vn- (B.l) 

The hamiltonian H is tri-diagonal in the basis (y? n )n>o and such that the /„ obey the 
recursion relation 

fn = —fn-2 - —fL-i , n>2, (B.2) 

On— 1 On— 1 

with /o = const and /i(t) freely specifiable. The positive numbers a n , n > 0, parameterize 
the hamiltonian, see Eqs. (13') and (17") in [25]. In the n = 2 equation in (B.2) we 
allowed ao > for later convenience; compared to ao = this modifies some aspects 
related to normalizability but does not affect the overall structure. The term 'supertask' 
derives from the fact that for suitable choices of a n and initial conditions the system 
can undergo spontaneous excitations with a pattern familiar from classical supertasks. 
That is, each excitation is aroused by a faster excitation of higher numbered excitations, 
ad infinitum, leading to severe violations of determinism. This can be achieved without 
sacrificing quantum mechanical principles, even in cases where the temporally conserved 
norm \\if)\\ is finite. The construction capitalizes on the fact that the differential form of the 
time evolution law as set by the Schrodinger equation id t ip = Hifj may under-determine the 
time evolution for systems with an infinite dimensional Hilbert space. It is only through 
the subtle injection of additional pieces of information in terms of boundary conditions and 
domains that the associated propagation kernel, t h- >■ e ltH , if well-defined, fully determines 
the time evolution. The countable infinity underlying (B.2) and the corresponding basis 
( ( / 9 n)n>o is essential in this context - any truncation of (B.2) at some n < N would lead to a 
closed system of differential equations which fully determines all functions /„, < n < N, 
including f x . It is therefore a 'supertask' to borrow information from infinity to specify 
f\{t) by some principle. 

The case where the a n grow without bound is referred to as an "accelerated supertask". 
In this situation the Schrodinger equation under-determines the time evolution and a wave 
function with normalizable initial condition will in general instantaneously evolve into a 
nonnormalizable one. For example with / n (0) = <5 nj o all f n (t) are nonzero for arbitrarily 
small t > 0, with 

f n (t) = -a Q a l ...a n _ l t n + 0(t n ~ 2 ), n>l. (B.3) 
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For a n growing faster than n this suggests that ^2 n>0 f n (t) 2 will diverge for arbitrarily 
small t > 0, even if (B.3) only provides the leading terms of an asymptotic expansion. If 
the f n (t) converge for t — > oo, the limiting values are dictated by (B.l) 

/o n even, 

aid?, . . . a n -i 

(B.4) 

0-10,3 ■ ■ ■ O n -2 j. , s , , 
/i(ooJ n odd. 

02^4 . . . 

In simple cases like ao = 0, a n = a, n > 1, normalizability of the wave function can be used 
as a criterion to restore uniqueness of the time evolution [25]. In the situation (B.3), (B.4) 
above normalizability is necessarily violated when a /i(oo) ^ 0. To see this note that 

:=^/nW 2 "2ao/o fdsh(s), (B.5) 

n>0 J ° 

is formally conserved on account of (B.2). The rearranging of the sums is legitimate only 
^ Sri>i o n f n (t)fn+i(t) converges absolutely for all t. In the t — > oo limit the latter sum 
contains denumerably many identical terms of the form /o a o/i(°°); so with / ^ the 
vanishing of a /i(oo) is a necessary condition for normalizability. In the following we are 
interested in situations where ao/i(oo) ^ and argue that, instead, Borel summability can 
be used as a criterion to restore uniqueness of the time evolution. To this end we make 
contact with the results of Section 3.2. 

Consider specifically the case 

a n = ^2n(2n + l)(2n + 2), n>l, (B.6) 

and set m = 1 and 

m p (t) = (4a) - P /4v^(^_l)n /p/2(v ^ t) ^ p = 2 mod 2 (B 7) 

Then the recursion relation (3.23) for the moments is mapped into (B.2) with the specific 
a n in (B.6) and f = 1. In particular (B.3) is mapped into m p {t) = j^yt p ^ 2 + 0(t p / 2 ~ 1 ), 
as required, and the asymptotics (B.4) is mapped onto (3.33). As seen in Section 3.2 
the formal power series in t obtained from (3.23) is asymptotic to the exact result for 
short times. It is also Borel summable and the Borel transform uniquely defines an rri2(t) 
for all times with a nonzero asymptotics 1712(00). Hence all moments m p (t) are uniquely 
determined and as shown in Section 3 they also coincide with the ones defined by the 
complex propagation kernel via (2.1). 
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Hence the Borel transform in this case augments precisely the piece of information oth- 
erwise supplied by the construction of the propagation kernel, the latter however cannot 
be achieved explicitly even in the simple case of a sextic anharmonic oscillator. The failed 
equivalence of (2.1) to the (z p )n t moments highlights that a differently constructed semi- 
group may define a different time evolution with a different t — > oo asymptotics. 

Although (B.6), (B.7) provides a mathematical isomorphism the conceptual interpretation 
of the m p and the f n is of course different. The m p are averages of the observables x p 
and are independent of the choice of basis in the underlying state space. Their dynamics 
(3.23) arises from instances of the observable flow equation (3.17) with respect to the 
stochastic time t. The dynamics of each fixed observable is driven by the quartic action 
via L = -e s / 2 He _s / 2 = d 2 x — d x Sd x , and only when a complete set of observables is 
considered is the effective dynamics of the set governed by the tri-diagonal hamiltonian 
underlying (B.2). The f n in (B.l), on the other hand, are the coefficients of the time 
dependent Schrodinger wave function with respect to a preferred basis, and as such do not 
directly qualify as observables. The time variable t refers to the physical time. Nevertheless 
the correspondence (B.7) cuts both ways. It could be used to design examples where the 
quantum mechanical supertasks can be solved by Borel transform. Conversely the under- 
determination in the Schrodinger dynamics (3.23) is ultimately the reason for the failure 
of the Parisi-Klauder conjecture in the present context. 
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